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Exchange interactions among defects in semiconductors are commonly treated within effective- 
mass theory using a scaled hydrogenic wave- function. However such a wave-function is only appli- 
cable to shallow impurities; here we present a simple but robust generalization to treat deep donors, 
in which we treat the long-range part of the wavefunction using the well established quantum de- 
fect theory, and include a model central-cell correction to fix the bound-state eigenvalue at the 
experimentally observed value. This allows us to compute the effect of binding energy on exchange 
interactions as a function of donor distance; this is a significant quantity given recent proposals to 
carry out quantum information processing using deep donors. As expected, exchange interactions 
are suppressed (or increased), compared to the hydrogenic case, by the greater localization (or de- 
localization) of the wavefunctions of deep donors (or 'super-shallow' donors with binding energy 
less then the hydrogenic value). The calculated results are compared with a simple scaling of the 
Heitler-London hydrogenic exchange; the scaled hydrogenic results give the correct order of mag- 
nitude but fail to reproduce quantitatively our calculations. We calculate the donor exchange in 
silicon including inter-valley interference terms for donor pairs along the {100} direction, and also 
show the influence of the donor type on the distribution of nearest-neighbour exchange constants at 
different concentrations. Our methods can be used to compute the exchange interactions between 
two donor electrons with arbitrary binding energy. 

PACS numbers: 03.67.Lx, 71.55.Cn 



I. INTRODUCTION 

Accurate estimations of exchange interactions among 
semiconductor defects are crucial in determining the 
magnetic properties of doped semiconductors near to a 
metal-insulator transition [l[ and in assessing the poten- 
tial of such defects for potential applications in quantum 
information processing @, However, the long-range 
tails of defect wavefunctions make it difficult or impossi- 
ble to treat these systems using fully ab initio techniques, 
so such calculations are generally performed within em- 
pirical models such as the well-established effective mass 
approximation 0, S @| ■ The quantum-mechanical prob- 
lem for a single defect is thereby reduced to that of a sin- 
gle electron moving in an effective medium determined 
by the effective mass and static permittivity of the host 
crystal; in the simplest form of the theory, the solution 
becomes that of a scaled hydrogen atom. The exchange 
between two such defects can then be obtained from 
well-established treatments of exchange in the hydrogen 
molecule, ranging from the simple Heitler-London model 
to the more sophisticated approach of Herring and 
Flicker [H, Q which exactly accounts for the two-electron 
correlations in the limit of large separations. 

However, a number of complications arise when consid- 
ering real defect systems. First, the minima of the con- 
duction band may not be isotropic; this is so in the impor- 
tant case of silicon, where there are six degenerate min- 



ima lying along the {100} directions in the Brillouin zone 
and each having a significant anisotropy (mi = 0.98m e , 
m 1 1 = 0.19?7i e ), where m e is electron mass. This is com- 
monly dealt with either by simply adopting an effective 
isotropic dispersion with an appropriate average effective 
mass (for example, m* = (mym^) 1 / 3 ), or by making a 
variational ansatz for the wavefunction of the anisotropic 
system, such as that of Luttinger and Kohn pj, [H, @ . 

Second, the existence of more than one minimum needs 
to be accounted for; this results in an effective-mass wave- 
function of the form 

VXr-) = ^«„F„(r)0 nfe - o (rO, (1) 

n 

where the sum runs over the different minima, F n is a 
slowly-varying envelope function, and <f> (r) = cxp(i/co • 
f)u r (r) is the Bloch function for the nth minimum. 

' nko v ' 

The coefficients a n arise from the coupling of the differ- 
ent band minima by short-range (large-wavevector) com- 
ponents of the potential; for an isolated substitutional 
donor in a perfect crystal, they correspond to an irre- 
ducible representation of the Td point group. In this 
paper we concentrate on the ground state of defects in 
silicon, for which the terms in equation |T]) correspond to 
the six {100} directions and the appropriate representa- 
tion is the identity representation: 

a (0) = -4(1,1,1,1,1,1)- (2) 
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The existence of these distinct parts to the wavefunc- 
tion gives rise to interference terms in properties of defect 
pairs as a function of their separation. Such terms were 
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originally considered in electron hopping and recombina- 
tion rates [Tol ITiT ] , and later included in exchange interac- 
tions fil.ll2j]. Other more recent calculations are discussed 
below. Note, however, that the expressions used in [l[ to 
average the exchange over the different minima assume 
a particular physical origin for the exchange (namely the 
direct electron-electron Coulomb interaction term) which 
does not, in fact, dominate the problem. Although its 
magnitude has approximately the same asymptotic scal- 
ing as the full exchange, by itself the direct term would 
lead to a ferromagnetic spin-spin interaction. The ob- 
served interaction is always antiferromagnetic, as it must 
be for a two-electron system [Hj]). 

Both the anisotropy and multiple-minima problems 
can be solved while remaining largely or entirely within 
the framework of effective-mass theory, and a number 
of recent papers have addressed the exchange between 
hydrogenic donors by this method. Much of the inspi- 
ration for this work came from the proposal ;2j to use 
such exchange in shallow donors as a qubit-qubit cou- 
pling mechanism to drive entangling gates in quantum 
information processing. The importance of precise donor 
positioning in determining the value of the exchange was 
emphasized in 13], but this paper follows [l| in using 
only the (ferromagnetic) direct exchange interaction for 
the calculations. The full Heitler-London formula for the 
exchange was used in a following paper [15j . where the 
contribution of strain to modifying the interaction was 
also discussed. These calculations were performed by ne- 
glecting the spatial variation of the periodic part of the 
Bloch functions u n ? appearing in equation fTJ) ; this ap- 
proximation was carefully examined in [l6| and found to 
be very accurate, and this finding was confirmed by a 
subsequent calculation of the electronic structure of de- 
fect pairs in which the Bloch functions from an ab initio 
treatment of the host crystal were used [17|. The same 
paper also found that allowing the phases of the different 
Bloch function contributions to 'float' relative to one an- 
other in the case of the defect pair made small difference 
to the results. 

The third difficulty is in many ways most serious: 
effective-mass theory predicts that all single-electron 
donors in a given host should have the same binding en- 
ergy, independent of their chemical nature. This is far 
from the case: in silicon, for example, observed binding 
energies range from 31.2 meV for Li and 45.6 me V for P 
through to 71.0 meV for Bi [Tq] , whereas the calculated 
value using the Luttinger-Kohn variant of effective mass 
theory is 31.3 me V It is very desirable to predict 

exchange in deeper donors by similarly simple methods: 
they are advantageous for quantum information process- 
ing Q because they are less prone to ionization, and 
have much longer spin-lattice relaxation lifetimes [1 91] — 
indeed, resonant Raman routes for spin relaxation are 
completely absent if the level spacing exceeds the maxi- 
mum phonon energy of the host material (64.5 meV for 
Si). 

Treating deep donors requires some significant correc- 



tions to effective mass theory. Two common approaches 
involve making an explicit short-range correction to the 
Coulomb field of the impurity (the so-called 'central cell 
correction'), or altering the long-range solution in the 
Coulomb field so that it corresponds to the observed 
binding energy (the 'quantum defect ' ap proach, which 
has its origins in atomic physics) [2(J, l2l| . However, the 
effect of these modifications on the exchange interactions 
of the defects has remained largely unknown. 

In this paper we show how the exchange in the most 
physically relevant range of separations can be com- 
puted between donors with arbitrary binding energies 
(including the important case of deep donors). We use 
a quantum-defect description, with a simple model po- 
tential to represent the central-cell correction; for rea- 
sons that are explained below, it is important to compute 
the exchange using a potential for which the wavefunc- 
tion concerned is an eigenfunction. We first introduce 
the techniques that we use: the effective-mass model, 
quantum-defect corrections to it and the simple central- 
cell corrections, followed by the methods we use to cal- 
culate exchange. Then we present our results for defects 
of different binding energies. 



II. METHOD OF CALCULATION 



A. Effective-mass theory 



The effective-mass equation 0, H, @] reads: 



[e n (k + -V) + U]F n = eF n 
i 



(3) 



where it is intended that the band energy e„ be expanded 
around the band extremum ko to second order-terms in 
(l/i)V. F n is the envelope function, in terms of which 
the true wavefunction is expanded using 



1> 



(4) 



where F n is a solution of the effective mass equation. 

In the simplest theory, the effective mass tensor is re- 
placed by a single averaged effective mass m* , resulting in 
an effective isotropic equation for the envelope function, 
which is then independent of the index n: 



2m* 



e 



e]F(r) = 0, 



(5) 



where e r is the relative permittivity of the host. In this 
paper we will follow other recent treatments (for exam- 
ple [H, [lj|) and work with this isotropic equation as our 
starting point. For silicon, m* — 0.33m e and e r = 11.7; 
this leads to a set of scaled atomic units for the hy- 
drogenic impurity problem (length = 1.94nm, energy 
Ha* = 0.062eV). 
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B. The quantum defect method 

The first component of our correction to effective-mass 
theory involves treating the wavefunction far from the 
impurity by the quantum defect method [21]. Even in 
the simplest isotropic approximation, the method was 
found to provide quantitative results giving good approx- 
imations to both the observed spectral dependence and 
magnitude of the photo- ionization cross section [2(| ■ An 
essential feature of the quantum defect method is that 
good approximate wave functions valid in the region out- 
side the impurity ion core can be determined using only 
a knowledge of the energy eigenvalues. Rather than at- 
tempting to solve it as an eigenvalue equation to deter- 
mine the allowed spectrum of e, equation §5§ is considered 
valid for large r only and solved for the asymptotic form 
of the envelope functions F(r) corresponding to the em- 
pirical value of e. Therefore the method can deduce the 
long-range part of the donor wave functions associated 
with the observed energy levels, whether deep or shallow, 
provided that the dominant corrections to effective-mass 
theory are short-range (operate only near the defect). Of 
the several different terms believed to contribute to the 
shift in binding energy for deep donors (see, for exam- 
ple, [12]), only the electron-phonon interaction operates 
far from the defect and its contribution is believed to be 
small. 

The most useful form of the far-field solution is a mul- 
tiple of the well-known Whittaker function which is just a 
particular linear combination of two standard linearly in- 
dependent confluent hypergeometric functions, the com- 
bination being determined by the boundary conditions at 
infinity. We consider the auxiliary radial function defined 
by P(r) = rR(r), and write 



P„,l(r) = N v ,,W Uii+ i(2r/u), 



(6) 



where W is a Whittaker function [23|], N is a normal- 
ization constant, and v is the quantum defect parameter 
defined by 



-1 
2^2 



(7) 



If and only if the donor is hydrogenic (i.e. shallow), v is 
an integer. 



C. Model central cell corrections 

For non-integer values of v, the radial function associ- 
ated with the solution ((6| diverges like R(r) = P(r)/r ~ 
- as r — > 0. (Nevertheless this divergence is integrable, 
so the functions have a well-defined normalization for all 
v [2l[.) In order to find a solution to equation (JSJ) that 
remains finite at the origin for a general given energy e, 
we have to correct the potential at short distances. One 
way to do this is to look for corrections based on the local 



physics of the impurity (for example incorporating cor- 
rectly the transition from a screened to an unscreened nu- 
clear potential 0, [2{| [26[ , or including self-consistently 
the scattering effects of the impurity by means of a pseu- 
doptential based on the microscopic physics [271]): an- 
other way is to correct the potential at small r empirically 
solely in order to make the solution regular there at the 
experimentally observed energy eigenvalue. In this sec- 
ond case the potential will not correspond to the physics 
operating in the core region of the real defect, but it will 
produce the correct shift in binding energy. We adopt 
the second approach in this paper, and refer to the em- 
pirical correction as a model central-cell correction. We 
should expect that the two approaches would give similar 
results for the calculation of exchange interactions, since 
they are determined predominantly by the long-range be- 
haviour of the wavefunction. 

We can then write the Hamiltonian of a single defect 
center A as 



A— iv- 



1 



Ra 



+ SV(r-R A ), 



(8) 



where the third term on the right hand side is model 
central cell correction, which is not unique. In this work 
we used two forms for SV: a <5-function shell or a square- 
well potential. 



1, S potential correction 

Let us choose SV to be a potential "shell" at radius 
r = a, where a is very small (a < 1 in scaled atomic 
units): 



V(r) 



1 



XS(r 



(9) 



In the following discussion, we will only consider s- 
symmetry wave functions because they are the donor 
ground states. After matching the large-r solution to a 
solution valid for small r that is regular at the origin and 
obeys the cusp condition for a nucleus of charge Z = 1 
(R'(0) = -R(0)), we find the solution is 



R(r) 



r < a 
r > a 



(10) 



where P u ,o{ r ) is the long-range solution for s-states de- 
fined by equation ([6J . Figure (JJ) shows the wave function 
for deep and shallow donors, both with the central cell 
correction and without it. Note how the diverging tails 
of the Whittaker function are cut off for r < a. 



2. Square-well potential correction 

An alternative matching is to replace the Coulomb po- 
tential by a square well, of depth Vq, for r < a: 



V(f) 



V 
i 



r < a 
r > a 



(11) 
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This model correction gives us: 
• E > V 



P(r) 

• E< V 

P(r) 



_ / i^,o(o) 



sin ka 



P»{r) 



r < a 
r > a 



r < a 
r > a 



(12) 



(13) 



In Figure $Q we compare the wave functions for both 
deep and shallow donors, with both 5-function square- 
well potential central cell corrections, and without any 
central-cell correction. Note that in both cases the di- 
verging tails of the wave function are cut off for r < a. 



D. Exchange calculations 

Within effective-mass theory, the exchange between 
two hydrogenic donors as a function of distance maps 
to the exchange between two hydrogen atoms as a func- 
tion of bond length. This problem has been treated using 
different methods since the 1920s. 



The Heitler- London model 




0.0 0.2 0.4 0.6 0.8 1.0 



(a) 




The Heitler-London model [7j simply evaluated the dif- 
ference between the expectation values of the Hamilto- 
nian in two two-electron wave functions, one a singlet 
( 1 S S ) and one a triplet ( 3 £ u ), both constructed from the 
Is ground states 4> a (r) and 4>b(r) of the single atoms. This 
approach neglects any contributions from other atomic 
states, and hence also contributions to the correlation 
arising from the polarization of one atom in the field of 
the other. The energy difference can be written as 



J = Ef — Ef, = 



0„(l)&(l)d[l]; 



J 
k 

■i 
3 

k' 



ft 2 - 2ks - k'], (14) 



(15) 
(16) 
(17) 
(18) 
(19) 



0„(1)[ ]0a(l)d[l]; 

rib 

&.(!)[-— ]&(l)d[l]; 
rib 

^(1)^(2)— d[l]d[2]; 
ri2 

0«(1)&(1)&(2)&(2)— d[l]d[2]. 

ri2 



All these quantities are positive; the one-electron contri- 
butions involving j and k yield a net antiferromagnetic 
coupling (positive J in our sign convention), while the 



(b) 

FIG. 1: Comparison of the wavefunctions generated by the 
S potential and square-well potential central-cell corrections, 
(a) Deep donor, v = 0.7; (b) super-shallow donor, v — 1.1. In 
each case, the dotted curve shows P„ (r)/r with a (^-potential 
central cell correction applied at a — 0.1, the dashed-dotted 
curve is P"^(r)/r with a square- well potential central cell 
correction applied with a = 0.1; dashed curve: uncorrected 
Whittaker function P^.o/r plotted for comparison; solid curve: 
radial part of a true Is Hydrogenic state (y = 1). Notice 
how the diverging tail of the Whittaker function is cut off for 
r < a by both central cell corrections; the form of the curves 
for r < a depends on the details of the potential applied. 



two-electron terms j' and k! are net ferromagnetic. How- 
ever, the neglect of correlation results in an unphysical 
logarithmic divergence in k! at large internuclear sep- 
arations R, causing the overall J to become negative 
(ferromagnetic), in violation of the theorem proved by 
Heisenberg 14| showing that the spatial ground-state of 
a two-electron system is always even under exchange of 
the particle positions, and hence is a spin singlet. 

We note that, in order to write J in the form of 
equation ([T4")) . it is essential that <p a and (j>b be ex- 
act eigenstates of the single- atom problem with energy 
Eq = —1/2. This allows one to replace the off-diagonal 
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matrix elements of the kinetic energy operator by 



0«(l)[-oV?]0 6 (l)d[l] 



&,(1)[£ + — ]<&(l)d[l]. 



(20) 



We refer to this as the 'Koiller method', since it is used eigenf unctions of the single-donor problem), an alterna- 

m mill, tive expression must be used: 

If (I20p is not obeyed (i.e. if cj> a and (pt, are not exact 



(V=f|£|iM = 



its 2 



[2(^a(l)|(-^V 2 )|0 a (l)) + 2(Mm-— )l^a(l)> 



1 



1 



+2<0 a (l)|( )|0 O (1)) T 2S«6,(1)|(--V')|&(1)> + 2<0 a (l)|( )|^(1))) 

Tib 2 ru 

d[1M2] mm t /• ^,(1^(1)^(2)^(2) i 

(</>-|i?|V->-(V>+W+>. 



(21) 
(22) 



We refer to this as the 'exact Heitler-London' method integral as 
below. 

In either case, when using a central-cell correction of 
the type discussed in i jll C[ it is important to include the 
central-cell terms in the appropriate matrix elements of 
the single-particle potential, i.e. in the calculation of the 
quantities j and k via equations (fit?)) and (jTTJ) , and in the 
one-electron terms in (l2~Tj) . 



2. Other trial wave] "unctions 

Other trial wavefunctions have been proposed for the 
hydrogenic case in an attempt to give some account of 
electron-electron correlation and remove the unphysical 
sign change in J. Kolos eq al. [29[ used trial wave func- 
tions in elliptic-coordinates similar to those proposed by 
[HI to perform a variational calculation of singlet and 
triplet energies, obtaining the exchange splitting from 
the energy difference. The numerical results are com- 
pared with the Heitler-London values in Figure ([2]). 



3. The Herring-Flicker asymptotic form 

A different approach was pursued by Herring [8] and 
by Gor'kov and Pitaevskii [3fJ , who pointed out that the 
exchange could be written in terms of a hyper-surface 



J=\( dS[(P$i)V$i - $iV(P$i)] +0(c 
2 Js 



-AR 



where 



), (23) 



(24) 



is a spatial wavefunction that is not properly antisym- 
metrized and has electron 1 localized on atom a and elec- 
tron 2 on atom 6, while (assuming the nuclear separation 
is in the z direction) the hypersurface S is defined by the 
condition Z\ = Z2- 

Herring and Flicker Q showed that this formula could 
be evaluated exactly in the limit R — > oo, giving an 
asymptotic form for the exchange of 



Jk-F — ^triplet — Psinglct 

= l.M2R b ' 2 e- 2R + 0{R 2 e- 2R ) 



(25) 



No unphysical sign change occurs in the Herring-Flicker 
asymptotic form. 

We compare the Heitler-London, Kolos and Herring- 
Flicker results for values of R up to 7ao (using a logarith- 
mic scale) in Figure [5Ja); we also give a comparison for 
larger values of R (up to 20ao) in Figure [U^b). Note that 
all three calculations agree quite well up to R = 7ao (the 
maximum radius for which Kolos' numerical results were 
calculated), and that the Heitler-London and Herring- 
Flicker forms are also close for R < 20ao. We therefore 
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expect that the Heitler-London approximation is a rea- 
sonable one in this distance range (the unphysical sign 
change to a ferromagnetic interaction does not occur un- 
til R = 49.5 ao). For comparison, Figure^b) also shows 
the variation of the quantity /ioff 2 / Airr 3 , which deter- 
mines the magnitude of the magnetic dipolar interaction; 
we see that the dipolar interaction starts to dominate 
over exchange for hydrogenic defects when R > 12a^. 



lnlJ/Ilu] 




FIG. 2: The single-envelope exchange splitting J for hy- 
drogenic donors calculated by three different methods as a 
function of inter-nuclei distance R and shown on a loga- 
rithmic scale, (a) Heitler-London exchange (solid curve), 
Herring-Flicker asymptotic (dotted curve) and Kolos' numer- 
ical calculation (dashed curve) for R < 7ao; (b) Heitler- 
London exchange (solid curve), Herring-Flicker asymptotic 
(dotted curve) and magnitude of magnetic dipolar interac- 
tion for R < 20ao. Notice that Heitler-London model agrees 
closely with the Herring-Flicker asymptotic form, and that 
from around R = 12ao the dipole-dipole interaction is larger 
than the exchange interaction. 



exchange formula (fl4|) should be replaced [l5| by 

K,K> 
2|„ |2 



X 



ot v \ J vti (R) cos(Av - kn) ■ R, 



where the pair of donors are at Ra — 0, Rb — R and 
R >> Oq (effective Bohr radius). The second sum (in 
square brackets) in equation (|26|) refers to the reciprocal- 
lattice expansion of the periodic Bloch functions, u v (f) = 

^2k cL fi etK r > an d K V} kf, are band minima points. The 
full expression for J„ M is given in the Appendix of [l5l |: 
in the isotropic effective-mass approximation where the 
envelope functions F n if) are the same for each minimum, 
and assuming that rapidly oscillating terms in the inte- 
grals (proportional to e*( fe< )_ ^ <j) ) r 5 where r is one of the 
integrated variables) are negligible, J MI/ can be replaced 
by the exchange J w computed using the radial functions 
derived from the Whittaker functions. The final expres- 
sion for the exchange is then 

j(r) = DE i4i 2 i4<i v( "' M ] ( 2? ) 

x\a v \ 2 \a^\ 2 J w {R) cos(fc„ - k^) ■ R. 

E. Fitting Whittaker function with Is Gaussian 

Once modified by our model central cell correction, 
Whittaker functions are valid solutions for the single- 
impurity problem for the whole range of r. To evalu- 
ate the integrals appearing in J W (R), we expand these 
solutions as a sum of Gaussians and use the analytical 
formulae in [3lj. We write 



RG(r) = J2A n G n (r,B n ) (28) 
G n (r,B n ) = e- B " r2 . 



In all the calculations presented here, we use ten Is 
Gaussian-type orbitals to fit each Whittaker function. 

In Figure ([3]) we show comparisons of such fits with the 
original wave functions. We see that we obtain a very 
good fit over the physically interesting distance range 
R < I2a* . 



4- Inter-valley effects 

The calculations presented so far give the exchange for 
molecular hydrogen. When we consider impurity states 
in a semiconductor we must remember that the full wave- 
function is of the form (j4|) . Assuming that the interval- 
ley couplings are mainly determined by single-impurity 
physics (so equation © still holds), the Heitler-London 



III. RESULTS 
A. Choice of parameters 

We select parameters appropriate for donors in silicon, 
as described previously: the appropriate effective mass 
is m* = (m^_ * mil) 3 = 0.33m e and the relative per- 
mittivity e r = 11.7, leading to an effective Bohr radius 
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P(r)/r 




(c) 

FIG. 3: Gaussian fits of Whittaker functions for (a) typical 
deep donor [y = 0.7), (b) typical hydrogenic shallow donor 
electron (y = 1.0) and (c) 'super-shallow' donor (y = 1.1). In 
each case the solid curve is the Gaussian fit and the dashed 
curve is the Whittaker function P(r)/r. 



a'X = 0.053nm e ' m = 1.94nm and an effective Hartree 

U m e 

Ha* = 27.2eV(a /a* Q ) = 62 meV. 

We take as a typical deep donor the case of Bi:Si, for 
which the binding energy is about 71meV [l8| or approx- 
imately 1.15 Ha*; according to formula ([7]), we therefore 
adopt v = 0.7 as a typical value for a deep donor, as 
used in Figure ([T]). For comparison and to show that our 
method is also capable of dealing with the case v > 1, we 
also show results for the case v = 1.1, corresponding to a 
hypothetical 'super-shallow' donor whose binding energy 
is less than the effective-mass value. (This is of less prac- 
tical interest, as the shallowest known donor, Li: Si, has 
v almost exactly equal to 1.) In all case where a model 
central-cell correction (of either type) is used, the radius 
is set to a = O.OlaQ. 

The results presented here cover the distance range 
from a*, to 12aj5- This range includes the lengthscales 



of greatest interest in practice, both for studies of the 
metal-insulator transition and for applications in quan- 
tum information processing, since it runs from the typi- 
cal nearest-neighbour separation at the highest attain- 
able densities below the metal-insulator transition (at 
a donor density tld — 4 x 10 18 cm -3 , the mean inter- 
donor spacing is 3.5 nm or 1.8 a*, and approximately 90% 
of nearest-neighbour separations are greater than a*,) up 
to separations where the dipole-dipole interaction starts 
to dominate over exchange (see Figure [2]). 



B. Benchmark: hydrogenic donors, v = 1.0 

First, we perform a benchmark calculation in which 
we compute the Heitler-London exchange interaction be- 
tween two pure hydrogenic donors (i.e. v — 1) by two 
different methods: first using the standard evaluation of 
the integrals from the exact Coulomb wavefunctions, and 
second using our Gaussian fit. The results are shown in 
Figure JU); the excellent agreement between the dotted 
line A (exact Coulomb wavefunctions) and the solid line 
C (Gaussian fit) shows that the errors introduced by the 
Gaussian fit are negligible. There is however some de- 
viation at large distances [R > 8<2q) between the solid 
curve C (exact Heitler-London result) and the dashed 
curve B (Koiller method); this arises because our Gaus- 
sian fit is not an exact eigenfunction of the one-centre 
problem and culminates in an unphysical cusp in curve 
B (corresponding to a sign-change in the exchange) at 
R « lla . 



C. Deep donors, v = 0.7 

Now we know that the errors arising from our Gaus- 
sian fits are likely to be small in the region of interest, we 
can move on to calculate our main result: the exchange 
coupling between two deep donors. The results for J w , 
calculated from a single envelope function, are shown in 
Figure ([5]). Notice that the six curves shown have very 
similar behaviours over most of the distance range. The 
three curves calculated using the 'exact Heitler-London 
method' and the different model central-cell corrections 
(A, C, E) are very close to one another, but in this case 
there are some significant deviations among the curves 
calculated using 'Koiller's mcthod'(B, D, F), both from 
the exact Heitler-London results and among one another. 
These deviations arise because the left side and right 
sides in equation (|20p are not exactly equal: in the cases 
of curves D and F, where central-cell corrections are in- 
cluded, this discrepancy arises from errors in the Gaus- 
sian fit, whereas in the case of curve B (no central-cell 
correction) there is an additional error because the Whit- 
taker function is not a true eigenfunction of the potential. 
This introduces a larger error and leads to an (unphysi- 
cal) cusp corresponding to a sign-change in the predicted 
exchange at R « Ha*, 
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FIG. 4: (Colour online.) The exchange coupling between 
hydrogenic donors as a function of inter-nuclear distance R. 
Blue dotted curve (A): exchange calculated directly from Hy- 
drogenic Is Coulomb wave function; red dashed curve (B): 
exchange from Koiller's method based on -R^oi black solid 
curve (C): exact Heitler-London method from Ri,q- 



Comparing the computed exchange with the results for 
hydrogenic states (Figured]), we see that the dominant ef- 
fect is the change in exponential decay constant reflecting 
the change in the exponent of the radial wavefunctions 
at large distances. It is therefore natural to ask whether 
a good approximation to exchange in deep donors can 
be obtained simply by re-scaling the results from the 
Heitler-London approach for a hydrogenic donor, mak- 
ing the replacement R — > Rjv. This comparison is made 
in Figure [UJa); we see that the scaled interactions have 
qualitatively the right behaviour and the correct order of 
magnitude, but do not match the details of our calcula- 
tion well. This impression is confirmed by looking at the 
ratio of the scaled to exact results shown in Figure Hfb): 
the error introduced by using the scaled approximation 
is nearly one order of magnitude over the distance range 
shown here, and increases still further at larger inter- 
donor distances. 



D. 'Super-shallow' donors, v = 1.1 

We also show the single-envelope exchange couplings 
J W (R) in the case of two 'super-shallow' donors. The 
results are shown in Figure 0. Once again the re- 
sults are independent of the type of central-cell correction 
used; furthermore (in the range R < 12ag) they are now 
largely independent of whether the exact Heitler-London 
approach or the 'Koiller method' is used . 



FIG. 5: (Colour online.) The exchange coupling between two 
deep donor electrons (y — 0.7) as a function of inter-donor 
distance R. Black solid curve (A): exact Heitler-London 
method without central cell correction; red dashed curve (B): 
Koiller's method without central cell correction; green dot- 
ted curve (C): exact Heitler-London method with 5-potential 
central cell correction; blue dashed-dotted curve (D): Koiller's 
method with 5-potential central cell correction; cyan dashed- 
dotted-dotted curve (E): exact Heitler-London method with 
square-well central cell correction; yellow short-dashed curve 
(F): Koiller's method with square- well central cell correction. 
(D and F are indistinguishable at long range). The cusp in 
curve B, corresponding to a sign change at R w llaS is un- 
physical (see text). 



Comparing the results to the hydrogenic case (Fig- 
ure S]), the dominant difference is once again the change 
in the exponent. As previously, we compare with the 
scaled hydrogenic Heitler-London result. Figure [8^a) 
shows that the matching is much closer than for the deep 
donors; the relative error is shown in Figure 0(b) . In 
contrast to the deep-donor case, the scaled hydrogenic 
wavefunction overestimates exchange interaction at short 
range, but underestimates it at long range. It is perhaps 
not surprising that the error is smaller in this case, since 
v is closer to the hydrogenic value v = 1 and so the scal- 
ing has a smaller effect. 



E. Inter-valley effects 

Having calculated J W (R), we include the inter- valley 
effects by using equation ([27]) . We have calculated the 
exchange between typical deep donors (y = 0.7), shallow 
donors (y = 1.0), and super-shallow donors (v = 1.1) 
in silicon in the case where the donor pair is along a 



FIG. 6: Comparison of exact Heitler-London exchange for 
deep donors (y — 0.7) with scaled hydrogenic results, (a) 
Exact Heitler-London exchange J OX act based on R§ 7 with a 
square-well potential central-cell correction (solid line) and 
scaled exchange splitting using the hydrogenic Heitler-London 
formula J sca ied = Jh-l , hydrogenic (7?/0.7) (dashed line), (b) 




The ratio 



as a function of inter-donor distance. 



R/a„ 



(100) direction from lap" to 12oq and both donors are 
in the ground state (A\ symmetry, so equation Q ap- 
plies) . We see the exchange interactions still decay expo- 
nentially over the whole range, but this decay is mixed 
with oscillations due to the inter-valley terms as shown 
in Figure ©. 

The solid points in the figure denote the distances cor- 
responding to the cubic lattice constant in Si (i.e., the ac- 
tual separations of substitutional sites along [100]). For 
all three defect types, the interference terms produce de- 
viations from site to site of approximately one order of 
magnitude, in agreement with the results for hydrogenic 
defects QJ]. 



FIG. 7: (Colour online.) The exchange coupling between 
two super-shallow donor electrons (u — 1.1) as a function of 
inter-donor distance R. Black solid curve (A): exact Heitler- 
London method without central cell correction; red dashed 
curve (B): Koiller's method without central cell correction; 
green dotted curve (C): exact Heitler-London method based 
on with 5-potential central cell correction; blue dashed-dotted 
curve (D): Koiller's method with 5-potential central cell cor- 
rection; cyan dashed-dotted-dotted curve (E): exact Heitler- 
London method with square-well central cell correction; yel- 
low short-dashed curve (F): Koiller's method with square- well 
central cell correction. Notice that these six curves have very 
similar behaviour over the whole range except at short dis- 
tance where there is a dependence on the central cell correc- 
tions. 



F. Distributions of exchange couplings 

In the interpretation of experiments on the ensem- 
ble of interacting spins, it is frequently the probability 
distribution of the nearest-neighbour exchange couplings 
which is the quantity of most relevance [l|, [32[ . In order 
to show how this depends on the donor type, we have 
used a continuous distribution [33[ of donor positions at 
two different densities, chosen so that the mean nearest- 
neighbour separations are respectively 3.5ag w 6.8nm 
and 7.0q,o ~ 13.5 nm, and plotted the distributions of 
log J for the three different types of defects in Figure (fT0|) . 
(The lower limit plotted corresponds approximately to 
the values of J where dipolar interactions begin to domi- 
nate over exchange.) All the distributions are extremely 
broad (as would be expected from the exponential scaling 
of the exchange with separation) but the strong depen- 
dence on the type of donor is evident. 

This distribution is important for the choice of donor 
concentration in samples for applications in quantum in- 
formation processing; for example, in the scheme pro- 
posed in [3j the density should ideally be sufficiently low 



that typical exchange interactions between neighbours 
produce a small evolution of the system on the timescales 
of gate operation. For deep donors with v — 0.7, typical 
nearest-neighbour interactions at the lower density are 
10 _4 Ha* or smaller, constraining gate operation times 
to be at most 0.1 ns. 



IV. CONCLUSION 

We have shown that a combination of quantum-defect 
theory (where Whittaker functions are used to describe 
simply the long-range part of donor wavefunctions) with 
model central-cell corrections can be used to describe the 
electronic structure of donors in semiconductors having 
binding energies significantly different from the ideal hy- 
drogenic effective-mass value, including in particular the 
important case of deep donors. From these wavefunc- 
tions we have shown that it is possible to calculate the 
exchange interactions between donors, both in the sin- 
gle envelope- function approximation and by including the 



10 




FIG. 8: Comparison of exact Heitler-London exchange for 
super-shallow donors (y = 1.1) with scaled hydrogenic results, 
(a) Exact Heitler-London exchange J CX act with a square-well 
potential central-cell correction (solid curve) and scaled ex- 
change splitting using the hydrogenic Heitler-London formula 
Jscaicd = Jh-l, hydrogenic (-R/L1) (dashed curve), (b) The ra- 
tio J ° calcd as a function of inter-donor distance. 



-5 - 



* -15 



C -20 



-25- 




R/a. 



interference effects between contributions from different 
conduction-band minima. These interference effects typ- 
ically cause the exchange to fluctuate by approximately 
one order of magnitude between the separations of suc- 
cessive substitutional sites, in agreement with previous 
calculations for shallow donors. 

From the comparison between the exchange calcula- 
tions performed with and without different central cell 
corrections, we can see that as expected the correct 
long-range behaviour of the wavefunction is more impor- 
tant than the form near the nucleus, except when the 
inter-donor distance becomes very small — however, the 
Heitler-London like approach we use for exchange is not 
expected to be accurate at short range. At long range, 
the central-cell corrections have only small effects on the 
exchange, although these may be enhanced by the in- 
terference of inter-valley terms. However, the inclusion 
of central-cell corrections is important to obtain consis- 
tent results for the exchange between the exact Heitler- 
London approach (equation (|2ip ) and the frequently- used 
simplification given by equation (|14p . 

Our calculations rely on a Gaussian fit to the true form 
of the radial function, which is least accurate very near 
the nucleus and at large distances. At large donor sep- 
arations we expect the error due to corrections near the 
nucleus to scale as 0(et 3 ), where a is the starting radius 
for the fitting. In our calculations a = O.Olag, much 
less than the donor separations 1 — 12a^ that we con- 
sider, so we expect this error to be negligible. The long- 
range error in the fit is important only at very large donor 
separations where the spin-spin interaction is no longer 
exchange- dominated. 

We have shown that the type of defect has a significant 
effect on the distribution of nearest-neighbour exchange 
couplings experienced at a given density. To summarize 
the magnitude of the effect, at at separation of 3.5 etQ « 
6.8 nm the exchange is approximately 0.15 meV for typ- 



FIG. 9: (Colour online.) Donor-pair exchange including inter- 
valley effects for pair separations along the {100} direction as 
a function of donor distance R. The exchange at integer mul- 
tiples of the lattice constant is represented by square points. 
A square-well central-cell correction and the exact Heitler- 
London approach were used for all the calculations. Lower 
(black) group (A, A'): two deep donors with v — 0.7; middle 
(red) group (B, B'): two typical shallow donors with v = 1.0; 
the upper (blue) group (C, C): two super-shallow donors 
with v — 1.1. 



ical deep donors, 3meV for hydrogenic shallow donors, 
and 5 meV for the so-called 'super-shallow' donors. This 
emphasizes the usefulness of deep donors for the short- 
term storage of quantum information, and confirms them 
as candidates for quantum information processing pro- 
vided that entangling interactions between them can be 
switched sufficiently quickly 0, [1[ . 

We have also shown that simply scaling the hydrogenic 
exchange interaction from the Heitler-London formula 
does not agree quantitatively with our explicit calcula- 
tions, though it does give the correct order of magnitude 
for the exchange. The error is most serious for the im- 
portant case of deep donors. 

We believe the methods we present here are quanti- 
tatively reliable, and simple enough to provide a useful 
tool for calculating exchange interactions between donor 
pairs having arbitrary binding energies. 
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FIG. 10: (Colour online.) The probability distribution of 
log J as a function of J (logarithmic scale) for donors with 
three different binding energies at densities 0.00396 (aJ5)~ 3 = 
5.42 x 10 17 cm -3 (curves A, B, C) and 0.0005 (aj) -3 = 
6.85 x 10 16 cm" 3 (curves D, E, F). Solid (black) curve A and 
dashed-dotted (blue) curve D: deep donor v — 0.7; dashed 
(red) curve B and dashed-dotted-dotted (magenta) curve E: 
shallow donor v = 1.0; dotted (green) curve C and short- 
dashed (yellow) curve F: super-shallow donor v = 1.1. 
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